/*==============================================================================
案例 6：房价预测分析（回归算法）
================================================================================

业务场景：
房地产公司需要建立房价预测模型，帮助定价决策和投资评估。
通过分析房屋特征（面积、位置、配置等）预测合理的房价。

学习目标：
1. 掌握非线性回归建模
2. 学习特征工程和变量转换
3. 理解正则化回归的应用
4. 掌握模型解释和业务应用

数据来源：auto.dta（模拟为房产数据）

作者：Stata ML Course
日期：2025-11-03
==============================================================================*/

clear all
set more off
capture log close

* 设置工作目录
cd "`c(pwd)'"

* 创建输出目录
capture mkdir "output/cases"
capture mkdir "output/cases/figures"
capture mkdir "data/cases"

* 开始日志记录
log using "output/cases/case06_house_price.log", replace text

display "=========================================="
display "案例 6：房价预测分析"
display "=========================================="
display ""

/*------------------------------------------------------------------------------
第一部分：数据准备（模拟房产数据）
------------------------------------------------------------------------------*/

display "第一部分：数据准备"
display "------------------"

* 加载数据
sysuse auto, clear

* 将汽车数据转换为房产数据的概念映射
* price -> 房价（万元）
* mpg -> 能源效率评分（转换为地段评分）
* weight -> 建筑面积（平方米）
* length -> 房龄（年）
* headroom -> 层高（米）
* trunk -> 车位数量
* foreign -> 是否学区房

rename price house_price
rename mpg location_score
rename weight building_area
rename length house_age
rename headroom floor_height
rename trunk parking_spaces
rename foreign school_district

* 调整数值范围使其更符合房产数据
replace house_price = house_price / 100  // 转换为万元
replace building_area = building_area / 10  // 转换为平方米
replace house_age = 30 - house_age / 5  // 转换为房龄（年）
replace location_score = location_score * 10  // 地段评分（0-100）

* 创建新特征
gen price_per_sqm = house_price * 10000 / building_area  // 单价（元/平米）
gen area_squared = building_area^2  // 面积平方项
gen age_squared = house_age^2  // 房龄平方项
gen luxury = (house_price > 60)  // 豪宅标记
gen new_house = (house_age < 5)  // 次新房
gen large_house = (building_area > 300)  // 大户型

* 创建区域变量（基于 rep78）
gen region = rep78
label define region_lbl 1 "郊区" 2 "远郊" 3 "中心区" 4 "核心区" 5 "CBD"
label values region region_lbl

display "样本量: " _N
display ""

* 数据概览
describe house_price building_area house_age location_score floor_height ///
    parking_spaces school_district

summarize house_price building_area house_age location_score

/*------------------------------------------------------------------------------
第二部分：探索性数据分析
------------------------------------------------------------------------------*/

display ""
display "第二部分：探索性数据分析"
display "------------------------"

* 1. 房价分布
histogram house_price, ///
    width(10) ///
    frequency ///
    title("房价分布") ///
    xtitle("房价（万元）") ///
    ytitle("频数") ///
    scheme(s2color)
graph export "output/cases/figures/case06_01_price_distribution.png", replace

* 2. 房价与建筑面积的关系
twoway (scatter house_price building_area, mcolor(blue%40)) ///
       (lfit house_price building_area, lcolor(red) lwidth(thick)) ///
       (qfit house_price building_area, lcolor(green) lwidth(thick)), ///
    title("房价与建筑面积的关系") ///
    xtitle("建筑面积（平方米）") ///
    ytitle("房价（万元）") ///
    legend(order(1 "实际数据" 2 "线性拟合" 3 "二次拟合")) ///
    scheme(s2color)
graph export "output/cases/figures/case06_02_price_vs_area.png", replace

* 3. 房价与房龄的关系
twoway (scatter house_price house_age, mcolor(blue%40)) ///
       (lfit house_price house_age, lcolor(red) lwidth(thick)), ///
    title("房价与房龄的关系") ///
    xtitle("房龄（年）") ///
    ytitle("房价（万元）") ///
    legend(order(1 "实际数据" 2 "线性拟合")) ///
    scheme(s2color)
graph export "output/cases/figures/case06_03_price_vs_age.png", replace

* 4. 不同区域的房价对比
graph box house_price, over(region) ///
    title("不同区域的房价分布") ///
    ytitle("房价（万元）") ///
    scheme(s2color)
graph export "output/cases/figures/case06_04_price_by_region.png", replace

* 5. 学区房 vs 非学区房
graph box house_price, over(school_district) ///
    title("学区房 vs 非学区房价格对比") ///
    ytitle("房价（万元）") ///
    scheme(s2color)
graph export "output/cases/figures/case06_05_school_district_effect.png", replace

* 6. 相关性矩阵
correlate house_price building_area house_age location_score floor_height parking_spaces

/*------------------------------------------------------------------------------
第三部分：数据分割
------------------------------------------------------------------------------*/

display ""
display "第三部分：数据分割"
display "------------------"

set seed 20251103
gen random = runiform()
gen train = (random <= 0.75)

count if train
display "训练集样本量: " r(N)
count if !train
display "测试集样本量: " r(N)

/*------------------------------------------------------------------------------
第四部分：模型 1 - 简单线性回归
------------------------------------------------------------------------------*/

display ""
display "第四部分：简单线性回归"
display "----------------------"

regress house_price building_area if train

estimates store model1

display ""
display "模型 1 结果："
display "每增加 1 平米，房价增加 " %6.2f _b[building_area] " 万元"
display "R² = " %5.3f e(r2)

/*------------------------------------------------------------------------------
第五部分：模型 2 - 多元线性回归
------------------------------------------------------------------------------*/

display ""
display "第五部分：多元线性回归"
display "----------------------"

regress house_price building_area house_age location_score ///
    floor_height parking_spaces school_district if train

estimates store model2

display ""
display "模型 2 结果："
display "R² = " %5.3f e(r2)
display "调整 R² = " %5.3f e(r2_a)

/*------------------------------------------------------------------------------
第六部分：模型 3 - 非线性回归（多项式）
------------------------------------------------------------------------------*/

display ""
display "第六部分：非线性回归"
display "--------------------"

regress house_price building_area area_squared house_age age_squared ///
    location_score floor_height parking_spaces school_district ///
    c.building_area#c.location_score if train

estimates store model3

display ""
display "模型 3 结果（包含非线性项）："
display "R² = " %5.3f e(r2)
display "调整 R² = " %5.3f e(r2_a)

/*------------------------------------------------------------------------------
第七部分：模型 4 - Lasso 回归
------------------------------------------------------------------------------*/

display ""
display "第七部分：Lasso 回归"
display "-------------------"

lasso linear house_price building_area area_squared house_age age_squared ///
    location_score floor_height parking_spaces school_district ///
    i.region c.building_area#c.location_score c.building_area#c.house_age ///
    c.location_score#i.school_district if train, ///
    selection(cv) rseed(20251103)

* 显示选择的变量
display ""
display "Lasso 选择的变量："
lassocoef, display(coef, standardized)

/*------------------------------------------------------------------------------
第八部分：模型对比与评估
------------------------------------------------------------------------------*/

display ""
display "第八部分：模型对比与评估"
display "------------------------"

* 在测试集上预测
foreach i in 1 2 3 {
    quietly estimates restore model`i'
    predict price_pred_m`i' if !train
}

* Lasso 预测
predict price_pred_lasso if !train

* 计算评估指标
display ""
display "模型性能对比（测试集）"
display "======================"
display ""
display "模型              RMSE      MAE       R²      MAPE"
display "---------------------------------------------------"

foreach model in m1 m2 m3 lasso {
    quietly gen sq_error_`model' = (house_price - price_pred_`model')^2 if !train
    quietly gen abs_error_`model' = abs(house_price - price_pred_`model') if !train
    quietly gen pct_error_`model' = abs((house_price - price_pred_`model') / house_price) if !train
    
    quietly summarize sq_error_`model'
    local rmse_`model' = sqrt(r(mean))
    
    quietly summarize abs_error_`model'
    local mae_`model' = r(mean)
    
    quietly correlate house_price price_pred_`model' if !train
    local r2_`model' = r(rho)^2
    
    quietly summarize pct_error_`model'
    local mape_`model' = r(mean) * 100
}

display "简单线性回归    " %7.2f `rmse_m1' "   " %7.2f `mae_m1' "   " %6.3f `r2_m1' "   " %6.2f `mape_m1' "%"
display "多元线性回归    " %7.2f `rmse_m2' "   " %7.2f `mae_m2' "   " %6.3f `r2_m2' "   " %6.2f `mape_m2' "%"
display "非线性回归      " %7.2f `rmse_m3' "   " %7.2f `mae_m3' "   " %6.3f `r2_m3' "   " %6.2f `mape_m3' "%"
display "Lasso 回归      " %7.2f `rmse_lasso' "   " %7.2f `mae_lasso' "   " %6.3f `r2_lasso' "   " %6.2f `mape_lasso' "%"
display ""

* 预测 vs 实际对比图（最佳模型）
twoway (scatter house_price price_pred_m3 if !train, mcolor(blue%40)) ///
       (function y=x, range(0 150) lcolor(red) lpattern(dash)), ///
    title("非线性模型：预测 vs 实际房价") ///
    subtitle("测试集，R² = " %5.3f `r2_m3') ///
    xtitle("预测房价（万元）") ///
    ytitle("实际房价（万元）") ///
    legend(order(1 "数据点" 2 "完美预测线")) ///
    scheme(s2color)
graph export "output/cases/figures/case06_06_prediction_vs_actual.png", replace

* 残差分析
quietly estimates restore model3
predict resid_m3 if !train, residuals

histogram resid_m3, ///
    normal ///
    title("残差分布（非线性模型）") ///
    xtitle("残差（万元）") ///
    ytitle("密度") ///
    scheme(s2color)
graph export "output/cases/figures/case06_07_residuals.png", replace

/*------------------------------------------------------------------------------
第九部分：业务应用
------------------------------------------------------------------------------*/

display ""
display "第九部分：业务应用"
display "------------------"

* 1. 房价预测示例
display ""
display "1. 房价预测示例"
display "----------------"

quietly estimates restore model3

* 场景 1：中心区，100平米，5年房龄，学区房
local pred1 = _b[_cons] + _b[building_area]*100 + _b[area_squared]*10000 ///
    + _b[house_age]*5 + _b[age_squared]*25 + _b[location_score]*70 ///
    + _b[floor_height]*3 + _b[parking_spaces]*1 + _b[school_district]*1 ///
    + _b[c.building_area#c.location_score]*100*70

display "场景 1 - 中心区学区房（100平米，5年，地段评分70）"
display "   预测房价: " %6.2f `pred1' " 万元"
display "   预测单价: " %8.0f (`pred1'*10000/100) " 元/平米"

* 场景 2：郊区，150平米，15年房龄，非学区房
local pred2 = _b[_cons] + _b[building_area]*150 + _b[area_squared]*22500 ///
    + _b[house_age]*15 + _b[age_squared]*225 + _b[location_score]*40 ///
    + _b[floor_height]*3 + _b[parking_spaces]*2 + _b[school_district]*0 ///
    + _b[c.building_area#c.location_score]*150*40

display ""
display "场景 2 - 郊区非学区房（150平米，15年，地段评分40）"
display "   预测房价: " %6.2f `pred2' " 万元"
display "   预测单价: " %8.0f (`pred2'*10000/150) " 元/平米"

* 2. 关键价格驱动因素
display ""
display "2. 关键价格驱动因素"
display "--------------------"

display "基于非线性模型的关键发现："
display "- 建筑面积：主要驱动因素，存在规模效应"
display "- 地段评分：每提高 10 分，房价增加约 " %5.2f (_b[location_score]*10) " 万元"
display "- 学区房溢价：约 " %5.2f _b[school_district] " 万元"
display "- 房龄折旧：每增加 1 年，房价降低约 " %5.2f abs(_b[house_age]) " 万元"

* 3. 投资建议
display ""
display "3. 投资建议"
display "------------"

* 识别被低估的房产（实际价格 < 预测价格）
gen price_gap = house_price - price_pred_m3 if !train
gen undervalued = (price_gap < -10) if !train  // 低估超过 10 万元
gen overvalued = (price_gap > 10) if !train  // 高估超过 10 万元

count if undervalued
display "潜在投资机会（被低估房产）: " r(N) " 套"
count if overvalued
display "高估房产（建议观望）: " r(N) " 套"

/*------------------------------------------------------------------------------
第十部分：保存结果
------------------------------------------------------------------------------*/

display ""
display "第十部分：保存结果"
display "------------------"

* 保存预测结果
preserve
keep if !train
keep house_price building_area house_age location_score school_district ///
    price_pred_m1 price_pred_m2 price_pred_m3 price_pred_lasso price_gap
export delimited using "output/cases/case06_predictions.csv", replace
display "预测结果已保存: output/cases/case06_predictions.csv"
restore

/*------------------------------------------------------------------------------
总结
------------------------------------------------------------------------------*/

display ""
display "=========================================="
display "案例 6 完成！"
display "=========================================="
display ""
display "主要成果："
display "1. 建立了 4 个房价预测模型"
display "2. 非线性模型表现最佳（R² = " %5.3f `r2_m3' "）"
display "3. 识别了关键价格驱动因素"
display "4. 提供了房价预测和投资建议"
display "5. 生成 7 个可视化图表"
display ""
display "关键发现："
display "- 建筑面积和地段是房价的主要驱动因素"
display "- 学区房有显著溢价"
display "- 房龄对价格有负面影响"
display "- 非线性模型能更好地捕捉价格规律"
display ""
display "输出文件："
display "- 图表: output/cases/figures/case06_*.png (7个)"
display "- 预测: output/cases/case06_predictions.csv"
display "- 日志: output/cases/case06_house_price.log"
display ""

log close

